
## Modified R script by Sofia Pagliarin 
## based on the syntax by the authors C. Kabengele & R. Hahn
## https://www.sciencedirect.com/science/article/pii/S0040162521003668
## Data available at doi:10.1016/j.techfore.2021.120934 

## Remove everything 
rm(list = ls())

## Set working directory
setwd("~/IHS/thesis/QCA-Analyse")
#dir()

## Import data file
 
 ## test import                                                                                
STREETEXPER <- read_delim("STREETEXPER-2.csv", 
delim = ";", escape_double = TRUE, locale = locale(date_format = "%d-%m-%Y", 
        decimal_mark = ",", grouping_mark = ".", 
        asciify = TRUE), na = "empty", trim_ws = TRUE)                                                                                                                                                               asciify = TRUE), na = "NA", trim_ws = TRUE)
View(STREETEXPER)
problems(STREETEXPER)

## Data management ***//*** ####

## Remove variables
names(STREETEXPER)
##STREETEXPER <- STREETEXPER[-(1),]

str(STREETEXPER)
#head(STREETEXPER)
names(STREETEXPER)
#View(STREETEXPER)
#STREETEXPER_qca$Cases # print first col


## List of variables

## Nr: no. 111
## Cases: instance street experiments in austrian and german cities 
## city : city where experiment happen
## EPX type of experiment 
## DIS % of space taken away from cars for other (transport ) uses
## DAYS : amount of days implemented after 16-03-2022 ( closing of austrian german border and uni's)
## PER [outcome] : State of experiment 1 it still there most recent observation on google earth
## AVG.2020 average grade score of AFDC farhradklima test 2020 (grade 6.0-1.0, 1 best grade)
##AVG.2022 AVerage grade score of AFDC farhradklima test 2022
##INF.FXX gradescore on specific questions of ADFC test (F1-F27) related to quality of bike infrastrure
##SAFE.FXX gradescore on specific questions of ADFC test (F1-F27) related to quality of bike safefy
##SAFE.WG worst score on safe failure condition
##ENVI-P/C pro-enviroment party (DIE GRÜNEN) largest party past(P)or current(C) with muncipal election 
##ENVI-SW year elections took place and new muncipal goverment formerd
##CAR : presence of car manufaction near city
## UNI.SK population of uni students in 1000 within city
## UNI.pop % of population in citie are students


## Descriptive statistics for each variable ***//*** ####
## Calibrating outcome not needed in cvs only 0 0.33 0.81 (temporary materials) 1.0(permant) are inserterd
# NA is unverivied data cases remove from data set
STREETEXPER.CLEAN <- STREETEXPER_2[complete.cases(STREETEXPER_2$PER, STREETEXPER_2$Name),]
##STREETEXPER.CLEAN <- STREETEXPER.CLEAN [-3,]
## variables with na are set to 0 not present
STREETEXPER.CLEAN$CAR[is.na(STREETEXPER.CLEAN$CAR)] <- 0
STREETEXPER.CLEAN$UNI.POP[is.na(STREETEXPER.CLEAN$UNI.POP)] <- 0
STREETEXPER.CLEAN$DIS[is.na(STREETEXPER.CLEAN$DIS)] <- 0
STREETEXPER.CLEAN$UNI.SK[is.na(STREETEXPER.CLEAN$UNI.SK)] <- 0
STREETEXPER.CLEAN$`ENVI-C`[is.na(STREETEXPER.CLEAN$`ENVI-C`)] <- 0
STREETEXPER.CLEAN$`ENVI-P`[is.na(STREETEXPER.CLEAN$`ENVI-P`)] <- 0
STREETEXPER.CLEAN$`ENVI-SW`[is.na(STREETEXPER.CLEAN$`ENVI-SW`)] <- 2019 # not data asumed no eceltions in or since 2020

summary(STREETEXPER.CLEAN)

sd(STREETEXPER.CLEAN$DIS)
sd(STREETEXPER.CLEAN$DAYS)
sd(STREETEXPER.CLEAN$AVG.2020)
sd(STREETEXPER.CLEAN$INF.AVG)
sd(STREETEXPER.CLEAN$SAFE.WG)
sd(STREETEXPER.CLEAN$CAR)
sd(STREETEXPER.CLEAN$PER)

##SUMMARY of data og code




## QCA 

## Calibration of variables into conditions ***//*** ####

##Load QCA package
#install.packages("QCA")
#install.packages("SetMethods")
library(QCA)
library(SetMethods)




##STREETEXPER$cPER<- round(calibrate(STREETEXPER$PER, thresholds = "e=.50, c=0.67, i=0.90" ), 2) 
# round , 2: rounded to 2 decimal points

## calbrarate infra grades 4,25 poor , 2,99 good (scale 6.0 to 1.0)
## make smart later
## try logistic = false to use lineair method
##INF.Grades <-  c(STREETEXPER.CLEAN, AVG.2020, INF.AVG, INF.1, INF.F3, INF.F18, INF.F20 ,INF.F22, INF.F25 )
##STREETEXPER.CLEAN$INFTEST <- round(calibrate(STREETEXPER[,INF.Grades], thresholds = "e=4.25, c=3.5, i=2.99" ), 2)
## CLACUlate and calubrate grade need to have 2 good grades (A*B)+(A*C)+(B*C)
##INF.K <- (STREETEXPER.CLEAN$INF.F1 + STREETEXPER.CLEAN$INF.F3)/2
##STREETEXPER.CLEAN$INF.K <- round(calibrate(INF.K, thresholds = "e=4.25, c=3.5, i=2.99" ), 2)
##INF.I <- (STREETEXPER.CLEAN$INF.F18 + STREETEXPER.CLEAN$INF.F20)/2
##STREETEXPER.CLEAN$INF.I <- round(calibrate(INF.I, thresholds = "e=4.99, c=4.0, i=3.25" ), 2)
##INF.N <- (STREETEXPER.CLEAN$INF.F22 + STREETEXPER.CLEAN$INF.F25)/2
##STREETEXPER.CLEAN$INF.N <- round(calibrate(INF.N, thresholds = "e=4.25, c=3.5, i=2.99" ), 2)
## minimal 2 good scores to be good infrastrure
##STREETEXPER.CLEAN$INF.Q<- pmax(pmin(INF.K,INF.I),pmin(INF.K,INF.N), pmin(INF.I,INF.N))
STREETEXPER.CLEAN$FsINF <- round(calibrate(STREETEXPER.CLEAN$INF.AVG , thresholds = "e=4.10, c=3.85, i=3.60" ), 2)
STREETEXPER.CLEAN$INF.AVG [which(STREETEXPER.CLEAN$INF.AVG <= 3.0)] <- 1.0
STREETEXPER.CLEAN$FsINF[which(STREETEXPER.CLEAN$FsINF == 0.5)] <- 0.5
## calibarate safety fail 5.0 and above pass is 4.0
STREETEXPER.CLEAN$SAFE <- round(calibrate(STREETEXPER.CLEAN$SAFE.WG , thresholds = "e=4.75, c=4.501, i=4.25" ), 2)
STREETEXPER.CLEAN$SAFE [which(STREETEXPER.CLEAN$SAFE == 0.5)] <- 0.51
STREETEXPER.CLEAN$SAFE [which(STREETEXPER.CLEAN$SAFE.WG >= 5.0)] <- 0.0


## calulation steps with NOR logic (A-B) Current - present
delta = STREETEXPER.CLEAN$`ENVI-C`- STREETEXPER.CLEAN$`ENVI-P`## delta == 0 continues goverment
year <- -2020 + (STREETEXPER.CLEAN$`ENVI-SW`) ## year 0 is 2020
ENV <- STREETEXPER.CLEAN$`ENVI-C` 
df <- data.frame(year,delta,ENV)
## no change since 2020   via formula year 0:4 
## if greens came to power sooner higher score
## ISSUE:  no idea how to implement a formula in R 
## 
## formula if green new to power 0.05 + 0.18(5-year) #switch 2020 -> 0.95
#df$ENV [which(df$delta == 1 & df$year >=0)] <- compute(0.05 + (0.18*(5 - (df$year))))]
## formula if greens out of power 0.95 - 0.18(5-year) -> 0.05
#df$ENV [which(df$delta == -1 & df$year >= 0)] <-  0.95 -(0.18*(5 - (df$year))) 

##df$ENV [(df$delta == 1 ) & (df$year >= 0)] <- 0.05 + (0.18*(4 - (df$year)))
##df$ENV [(df$delta == -1 ) & (df$year >= 0)] <-  0.95 -(0.18*(4 - (df$year))) 

##6 step method 
df$ENV [which(df$delta == 1 & df$year == 0)] <- 0.9
df$ENV [which(df$delta == 1 & df$year == 1)] <- 0.81
df$ENV [which(df$delta == 1 & df$year == 2)] <- 0.67
df$ENV [which(df$delta == 1 & df$year == 3)] <- 0.33
df$ENV [which(df$delta == 1 & df$year == 4)] <- 0.18
df$ENV [which(df$delta == -1 & df$year == 0)] <- 0.18
df$ENV [which(df$delta == -1 & df$year == 1)] <- 0.33
df$ENV [which(df$delta == -1 & df$year == 2)] <- 0.67
df$ENV [which(df$delta == -1 & df$year == 3)] <- 0.81
df$ENV [which(df$delta == -1 & df$year == 4)] <- 0.9
df$ENV [which(df$delta == 0 | df$year <= -1)] <- STREETEXPER.CLEAN$`ENVI-C`
STREETEXPER.CLEAN$ENVI <- df$ENV
STREETEXPER.CLEAN$ENVI[which (STREETEXPER.CLEAN$ENVI > 1)] <-1
STREETEXPER.CLEAN$ENVI[which (STREETEXPER.CLEAN$ENVI == 0.5)] <-0.67
STREETEXPER.CLEAN$ENVI[which (STREETEXPER.CLEAN$ENVI < 0)] <-0

## calibrate UNI 6 step
## issues with end condition
## Calibrate uni-pop as UNITOWN
##STREETEXPER.CLEAN$UNITOWN <- round(calibrate(STREETEXPER.CLEAN$UNI.POP , thresholds = "e=10, c=15, i=20" ), 2)
STREETEXPER.CLEAN$UNI <- 0
STREETEXPER.CLEAN$UNI [which(STREETEXPER.CLEAN$UNI.SK < 15 )] <- 0
STREETEXPER.CLEAN$UNI [which(STREETEXPER.CLEAN$UNI.SK >= 15)]  <- 0.33
STREETEXPER.CLEAN$UNI [which(STREETEXPER.CLEAN$UNI.SK >= 25 | STREETEXPER.CLEAN$UNI.POP >= 12)] <- 0.67
STREETEXPER.CLEAN$UNI [which(STREETEXPER.CLEAN$UNI.SK >= 30 & STREETEXPER.CLEAN$UNI.POP >= 16)] <- 1


##calibrate dis space
STREETEXPER.CLEAN$SPACE <- round(calibrate(STREETEXPER.CLEAN$DIS , thresholds = "e=15, c=20, i=49" ), 2)
STREETEXPER.CLEAN$SPACE[which(STREETEXPER.CLEAN$SPACE == 0.5)] <- 0.51
## reaction time quicker then 180 days
STREETEXPER.CLEAN$DAYS[is.na(STREETEXPER.CLEAN$DAYS)] <- 365
STREETEXPER.CLEAN$REACT <- round(calibrate(STREETEXPER.CLEAN$DAYS , thresholds = "e=180, c=60, i=30", logistic = TRUE), 2)
STREETEXPER.CLEAN$REACT [which(STREETEXPER.CLEAN$REACT == 0.5)]<- 0.51

##BAHN
##STREETEXPER.CLEAN$BAHN[which(STREETEXPER.CLEAN$BAHN == 1)] <- 0.33
STREETEXPER.CLEAN$BAHN[which(STREETEXPER.CLEAN$BAHN == 2)] <- 0.67
STREETEXPER.CLEAN$BAHN[which(STREETEXPER.CLEAN$BAHN >= 3)] <- 1


names(STREETEXPER.CLEAN)
Miss <- is.na(STREETEXPER.CLEAN[,c(9,28,31:37)])
print(Miss)
## remove imcomplete cases with na in conditions
STREETEXPER.QCA <- STREETEXPER.CLEAN[complete.cases(STREETEXPER.CLEAN,c(9,28,31:37))]
## Save QCA dataset 
names(STREETEXPER.QCA)
##STREETEXPER_qca <- STREETEXPER.CLEAN[ , c(24, 27:32, 5)] # outcome placed as last column
          ## 7 conditions, 1 outcome condition, 110 cases 


## Analysis of necessity ***//*** ####
## Performed for each individual condition
## Goal: checking superset relations 

##cnvector <- as.numeric(cnvector)
names(STREETEXPER.qca)
View(STREETEXPER.qca)
write.csv(STREETEXPER.QCA, "QCAFUZZYDATA.csv")
STREETEXPER.QCA <- QCAFUZZYDATA
View(STREETEXPER.QCA)
STREETEXPER.QCA$SAFE[which(STREETEXPER.QCA$SAFE == 0.50)] <-0.49
STREETEXPER.QCA <- STREETEXPER.QCA[,-c(1,4)]
conds <- subset(STREETEXPER.QCA, select=c(2:8)) # select conditions
conds <- as.numeric(conds)
names(conds)

nec_y <- QCAfit(conds, STREETEXPER.QCA$PER, names(conds), necessity = TRUE)

nec_y
#nec_y <- QCAfit(conds, 'PER',  relation ="nec" ) # SetMethods package
##nec_y <- pof(conds, 'PER',STREETEXPER_qca, relation = "nec") # QCA package
#nec_y # SIZE is necessary for OUT=1: consistency necessity > 0.9, high RoN, meaningful also theoretically 
##pof("ENVI <= PER",data = STREETEXPER_qca)

## Alternatively, it is possible to perform the analysis of necessity for each 
# individual conditions separately:
#pof("MNO <= PEN", data = STREETEXPER_qca) # pof: parameters of fit
#pof("~MNO <= PEN", data = STREETEXPER_qca) # tilde: ~ indicates absence of condition
#pof("GROUP <= PEN", data = STREETEXPER_qca)
#pof("~GROUP <= PEN", data = STREETEXPER_qca)
#pof("ATMDENS <= PEN", data = STREETEXPER_qca)
#pof("~ATMDENS <= PEN", data = STREETEXPER_qca)
#pof("REG <= PEN", data = STREETEXPER_qca)
#pof("~REG <= PEN", data = STREETEXPER_qca)
#pof("MOB <= PEN", data = STREETEXPER_qca)
#pof("~MOB <= PEN", data = STREETEXPER_qca)
#pof("AGENT <= PEN", data = STREETEXPER_qca)
#pof("~AGENT <= PEN", data = STREETEXPER_qca)
#pof("USE <= PEN", data = STREETEXPER_qca)
#pof("~USE <= PEN", data = STREETEXPER_qca)

## Repeat the performance of the analysis of necessity for the 
## non-occurrence (absence) of the outcome (i.e., not high penetration rate of mobile money)
nec_ny <- QCAfit(conds, STREETEXPER.QCA$PER, names(conds), 
                 necessity = TRUE, neg.out = T)
nec_ny

# nesity combination of conditions

# ISSUE ; nessity relantion doenst work coerced type double 
## Dusa h 5.5
QCAFUZZYDATA <- as.numeric(STREETEXPER.QCA)
## doesnt work double data type cannot be converted
##superSubset(QCAFUZZYDATA, outcome = PER, conditions = "FsINF, SAFE,UNI, ENVI, SPACE, REACT" ,incl.cut = 0.85)
##superSubset(QCAFUZZYDATA, outcome = ~PER,  ,incl.cut = 0.85)


## Analysis of sufficiency ***//*** ####
## Goal: check subset relations 

pof(ENVI, PER, data=QCAFUZZYDATA, relation = "suf")
pof(~ENVI, ~PER, data=QCAFUZZYDATA, relation = "suf")
##pof(SAFE, ~PER, data=QCAFUZZYDATA, relation = "suf")
##pof(~SAFE, ~PER, data=QCAFUZZYDATA, relation = "suf")
pof(~CAR, ~PER, data=QCAFUZZYDATA, relation = "suf")
pof(~INF.CQ, PER, data=QCAFUZZYDATA, relation = "suf")
pof(~REACT, PER, data=QCAFUZZYDATA, relation = "suf")
#tutoial dusa
#data("LF")
# using (LF, sum(fuzzyand(DEV,SURV))/sum(DEV))
pof("SAFE + ~SAFE", PER, data=QCAFUZZYDATA, relation = "suf")
pof("SAFE + ~SAFE", ~PER, data=QCAFUZZYDATA, relation = "suf")
pofind(data = QCAFUZZYDATA,outcome = "PER",conditions = "CAR,SPACE, REACT, SAFE,ENVI,UNI", relation = "nec")
##pof("REACT*~ENVI*~SAFE*UNI*SPACE",PER,data=QCAFUZZYDATA, relation = "suf")
pof("~REACT*~ENVI*~SAFE*UNI*SPACE + SPACE*REACT*~ENVI*~SAFE*~UNI",~PER,data=QCAFUZZYDATA, relation = "suf")pof(SPACE)
## Assign rownames as a combination of case name with outcome set membership value
city.names <- substr(STREETEXPER.CLEAN[[2]],start = 1, stop = 4)
exper.names <-substr(STREETEXPER.CLEAN[[4]],start = 1, stop = 5)

##QCAFUZZYDATA <- QCAFUZZYDATA[ , -c(3)]
QCAFUZZYDATA <- cbind(STREETEXPER.CLEAN$city,STREETEXPER.CLEAN$EPX,STREETEXPER.CLEAN$`area - street`,QCAFUZZYDATA)
write.csv(QCAFUZZYDATA, "STREETEXPERRAW.csv")
##rownames(STREETEXPER_qca) <- paste(rownames(STREETEXPER_qca),STREETEXPER.CLEAN[[5]],city.names,exper.names,sep="-")
##rownames(STREETEXPER_qca) <-paste(rownames,city.names, sep = "|")
##rownames(QCAFUZZYDATA) <- paste( exper.names, city.names, sep = "-")
QCAFUZZYDATA$SPACE[which(QCAFUZZYDATA$SPACE == 0.5)] <--0.51
## Generate truthtable for the occurrence of the outcome (Y)
tt_STREETEXPER <- truthTable(STREETEXPER.QCA, outcome = "PER", 
                      conditions = "CAR,SPACE, REACT, SAFE,ENVI,UNI,FsINF",
                      incl.cut = 0.81, complete = FALSE, show.cases = TRUE, 
                      sort.by = "n", n.cut = 3)
tt_STREETEXPER



## Conservative solution for the occurrence of the outcome
c_soly <- minimize(tt_STREETEXPER, details = TRUE)
c_soly

## Correspondence of the pathways ('causal recipes') to Table 2 in the article: 
## S2: REG*AGENT*~ATMDENS*MNO*GROUP +
## S3: MOB*AGENT*~ATMDENS*MNO*GROUP +
## S1: AGENT*~ATMDENS*MNO*GROUP*~USE +
## S4: REG*MOB*AGENT*ATMDENS*~MNO*GROUP*USE
## -> PEN

## A way to factorise the 4 solution pathways
## As per Fiss Chart in Table 2:
## MNO*GROUP*AGENT(~ATMDENS*~USE + REG*~ATMDENS + MOB*~ATMDENS) +
## ~MNO*REG*MOB*AGENT*ATMDENS*GROUP*USE
## -> PEN

STREETEXPER.QCA <- QCAFUZZYDATA[complete.cases(QCAFUZZYDATA,c(1:9))]
QCAFUZZYDATA[which(QCAFUZZYDATA$SAFE == 0.50)] <- 0.55
## Generate truthtable for the non-occurrence of the outcome (y or ~Y)
tt_STREETEXPER_ny <- truthTable(QCAFUZZYDATA, outcome = "~PER", 
                      conditions = "SPACE,FsINF,REACT,ENVI,CAR",
                      incl.cut = 0.75, complete = FALSE, show.cases = TRUE, 
                      sort.by = "n, incl", n.cut = 3, dcc = TRUE)
tt_STREETEXPER_ny


##QCAFUZZYDATA [,3] <- NULL
QCAFUZZYDATA$SAFE[which(QCAFUZZYDATA$SAFE == 0.50)] <- 0.55
## Conservative solution for the non-occurrence of the outcome
tt_STREETEXPER_ny <- truthTable(QCAFUZZYDATA, outcome = ~PER,conditions = "ENVI,REACT,SAFE,SPACE", incl.cut = 0.65, n.cut = 2)
tt_STREETEXPER_ny

## hypotesitng
c_solny <- minimize(tt_STREETEXPER_ny, details = TRUE, )
c_solny
## U1: ~REG*~MOB*~AGENT*~ATMDENS*~MNO*~GROUP +
## U2: ~REG*~MOB*~AGENT*~ATMDENS*GROUP*~USE + 
## U7: REG*MOB*~AGENT*~ATMDENS*~MNO*~GROUP + 
## U6: REG*MOB*~AGENT*MNO*GROUP*~USE +
## U3: REG*~AGENT*~ATMDENS*~MNO*~GROUP*~USE + 
## U4: ~MOB*~AGENT*~ATMDENS*~MNO*~GROUP*~USE + # From M2
## U5: ~REG*~MOB*~AGENT*~ATMDENS*GROUP*~USE + # From M2
## -> ~PEN 



## The performance of QCA implies also other steps, which we have not seen here. 
## Examples of such components of QCA are: 
## the solving of true logical contradictions, the generation of the 
## parsimonious and intermediate solutions (besides the conservative solution we did in this exercise)
## for the occurrence and non-occurrence of the outcome, 
## the definition of directional expectations, the evaluation of simplifying assumptions (easy and difficult counterfactuals),
## the generation of xy plots of the solution pathways and overall solution(s)
## and their visual inspection to "go back to the cases"
## (and in case to perform post-QCA selection of cases).

## This exercise was hence only a "preview" of the workings and goals of QCA by replicating 
## a QCA published study. 

## Full performance of QCA would require "going full circle" from sample construction to 
## the interpretation of the cases based on the obtained solutions. 




